Domino plasmons for subwavelength terahertz circuitry 
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A new approach for the spatial and temporal modulation of electromagnetic fields at terahertz fre- 
quencies is presented. The waveguiding elements are based on plasmonic and metamaterial notions 
and consist of an easy-to- manufacture periodic chain of metallic box-shaped elements protruding out 
of a metallic surface. It is shown that the dispersion relation of the corresponding electromagnetic 
modes is rather insensitive to the waveguide width, preserving tight confinement and reasonable ab- 
sorption loss even when the waveguide transverse dimensions are well in the subwavelength regime. 
This property enables the simple implementation of key devices, such as tapers and power dividers. 
Additionally, directional couplers, waveguide bends, and ring resonators are characterized, demon- 
strating the flexibility of the proposed concept and the prospects for terahertz applications requiring 
high integration density. 



I. INTRODUCTION 

Electromagnetic radiation in the terahertz (THz) 
regime is a central resource for many scientific fields such 
as spectroscopy, sensing, imaging, and communications. 
We are curr ently witnessing the take off of THz tech- 
nologies^! 3 -' with applications as diverse as astronomy^, 
medicine^, or security^. Within this general endeavour, 
the building of compact THz circuits would stand out as 
an important accomplishment. This requires the design 
of THz waveguides carrying tightly confined electromag- 
netic (EM) modes, preferably with subwavelength trans- 
verse cross section. Besides circuit integration and com- 
pact device design, sub-A localization may be advanta- 
geous for waveguide THz time-domain spectroscopy^ and 
non-diffraction-limited imaging^. Although a number of 
structures has been put forward, none of them passes all 
the following requirements. First, structures should be 
easily manufactured and, if possible, planar and mono- 
lithic. Second, as above motivated, subwavelength trans- 
verse size is also needed. Finally, absorption and bend 
losses should be small, and the waveguides sufficiently 
versatile for the design of key functional devices. Partic- 
ularly important are in/out-couplers since they work as 
the interface to external waves. In this context, compact 
tapers able to laterally compress the modes down to the 
sub-A level seem essential. 

With these requirements in mind it is illuminating to 
compare optical waveguides with those proposed for the 
THz range. Optical fiber modes exhibit a cut-off diame- 
ter such that fibers cannot be of subwavelength size and, 
moreover, modal size grows as this cut-off is approached. 
In the THz regime dielectric wires^meet the same limita- 
tion. The photonic crystal fiber approach^ shares the in- 
ability of sub-A confinement, and low-index discontinuity 
waveguides^ perform better regarding modal size but are 
neither planar nor monolithic. Metals are an alternative 



when operating at optical frequencies because surface 
plasmon polaritons (SPPs) propagating on flat metal- 
dielectric interfaces provide subwavelength confinement 
in the direction perpendicular to the surface. This metal- 
lic route has been tried for the THz regime, but the corre- 
sponding EM modes, known as Zenneck or Sommerfeld 

wa ves, l oose their confined character at these frequen- 
ciesl™ 

Geometrically-induced surface plasmons 1 supported 
by periodically corrugated metallic surfaces overcome the 
low localization of Zenneck waves. A number of waveg- 
uiding st ructures b ased on this concept has been already 
suggeste d 15 * 16 ' 17 * ^1, but they have failed to meet simul- 
taneously all the above mentioned requirements. In this 
work we present structures based on this metamaterial 
approach that, thank to their superior properties, ful- 
fill the requisite easy fabrication, subwavelength confine- 
ment, low loss, and device flexibility. 



II. RESULTS 

The basic structure consists of a periodic arrangement 
of metallic parallelepipeds standing on top of a metallic 
surface and resembling a chain of domino pieces, see in- 
set of Fig. [lja) . In contrast with corrugated wires^, this 
is a planar and monolithic system and should not pose 
significant manufacturing problems, its geometry being 
much simpler than corrugated V-groove*^ or wedges^. 
The properties of its guided modes, hereafter referred to 
as domino plasmons (DPs), are mainly controlled by the 
geometric parameters defining the dominoes, i.e., peri- 
odicity d, parallelepiped height h, lateral width L, and 
inter-domino spacing a. 

First, the basic properties of DPs are described starting 
with the dispersion relation. To gain physical insight we 
model the metal as a perfect electric conductor (PEC). 
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FIG. 1: Modal properties of domino plasmons. a, Dis- 
persion relation of DPs for various lateral widths L. Black and 
grey (blue) lines correspond to height h = 1.5d (h = 0.75d). 
Dashed line stands for infinitely wide dominoes (L = oo). 
Inset: diagram of the domino structure and geometric pa- 
rameters (the arrow depicts the mode propagation direction), 
b and c, Modal shape of DPs: transverse (xy) electric field 
(arrows) and horizontal (xz) electric field (color shading) for 
DPs of height h = 1.5d and widths L = 0.5d (b) and L = 8d 
(c). d and e, Same as in panels b and c, but now for a ID 
array of free-standing metallic rods (h = 3d). The designa- 
tions even and odd label the symmetries of the modes with 
respect to the corresponding white lines. The fields in b-d arc 
plotted for d/X = 0.125, marked with a red open dot in a, the 
white bar in c being the wavelength (valid for panels b-d). 
The field in e is computed for the same k as in panels b-d, 
corresponding now to d/X = 0.056. In panels a-e metals are 
modelled as PECs. f, DP modal effective index as a function 
of lateral dimension L in units of wavelength. Various oper- 
ating frequency regimes are considered: A = 1.6mm (red), 
A = 0.16mm (green), A = 0.016mm (blue), and A = 1.5/xm 
(magenta). To compute panel f, a realistic description of the 
metals is used. As described in the main text, the period- 
icity d is different for the various operating frequencies, and 
h = 1.5c/, a = 0.5d, L = 0.5d, 24d. 



Later, realistic metals will be considered, letting us to 
discern what is the frequency regime where DPs display 
the most promising behaviour. The numerical technique 
is described in the Methods section. Due to the simple 
scaling properties of PECs, we choose the periodicity d 
as the length unit. The value of a is not critical for the 
properties of DPs and is set as a = 0.5c? in this work. DP 
bands present a typical plasmonic character, i. e., they ap- 
proach the light line for low frequencies and reach a hor- 
izontal frequency limit at the edge (fc c d g0 = 7r/d) of the 
first Brillouin zone (Fig.[IJa), notice that only fundamen- 
tal modes are plotted). In all cases DP modes lie outside 
the light cone, a fact that explains their non-radiative 
character. While the limit frequency of SPPs for large 
k is related to the plasma frequency, the corresponding 
value for DPs is controlled by the geometry. The influ- 
ence of the height h is clear: the band frequency rises for 
short dominoes (h — 0.75c?, blue line) as compared with 
that of taller ones (h — 1.5d, black line). This can be un- 
derstood as follows. When L — oo, dominoes become a 
one-dimensional (ID) array of grooves which supports a 
geometrically-induced plasmon mode. Its dispersion re- 
lation is represented, for h = 1.5c?, with a dashed line in 
Fig. [TJa). An approximation for the dispersion relation 
in this limit, neglecting diffraction effects and for A >> d. 



k = ko\ 1 



tan 2 (f/j,/i), 



(1) 



where k is the modal wave vector, k = 2ir/\, and the y 
component of the wave vector inside the grooves, q y , is 
related with the corresponding x and z components by 
Qx + Qy + 1z = ■ The EM fields inside the grooves are 
independent of x and z for the chosen longitudinal polar- 
ization and thus q x = q z = 0, q y = ko- This means that 
Eq. [T] presents an asymptote for k^h = 7r/2 or, in other 
words, for a normalized frequency d/Ah. Such depen- 
dence with h explains the behaviour observed in Fig. [lja) . 
For h — 1.5c?, the value chosen in the rest of the paper, 
this estimation for the limit normalized frequency pro- 
vides a figure of about 0.17, not far from the numerically 
computed value of 0.14. 

The most striking characteristic of DPs is their be- 
haviour when the lateral width L is changed. All bands 
in the range L — 0.5c?, . . . , 3c? lie almost on top of each 
other (Fig.[TJa), grey curves). In other words, the modal 
effective index, n c ff = k/ko, is rather insensitive to lat- 
eral width. If we now look at the numerically computed 
dispersion relations of Fig. [l|a) in the light of Eq. [l] it 
seems that q y w fco for a wide range of L. Remarkably, 
the bands remain almost unchanged even for L = 0.5c?, 
whose modal size is well inside the subwavclength regime, 
as observed in the field distribution plotted in panel (b) 
for a narrow structure (L — 0.5c?). Notice that the white 
bar representing A in panel (c) is valid for panels (b)- 
(d), and the modal size, defined as the FWHM of the 
modulus of the Poynting field distribution, is 0.21 A for 
panel (b). The described behaviour is to be contrasted 
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with that of conventional plasmonic modes in the optical 
regime for which sub- A lateral confinement is not a trivial 
issued When the lateral size of the structure support- 
ing such modes becomes subwave length, either the modal 
size grows, or the effective index increases. This is best il- 
lustrated for a free standing metallic stripe 2 *. A stripe of 
large lateral width essentially supports an SPP mode but, 
as this width is decreased and becomes sub- A, the plas- 
monic modes that are not cut-off become either highly 
confined but quite lossy, or otherwise highly delocalized. 

The insensitivity of DPs to lateral size, which consti- 
tutes the foundation of some of the devices presented 
later, can be linked in an intuitive way to the structure 
geometry. The key role played by the metallic substrate 
can be explained by comparison of DPs and the modes of 
a Yagi-Uda struc ture, i.e ., a periodic array of free stand- 
ing metallic rod d 22 * 23 * 2 ^. For small L the dispersion of 
DPs of height h is identical to that of the fundamen- 
tal mode of a Yagi-Uda structure of height 2h, and the 
field of the DP is simply the upper half of the funda- 
mental mode of a rod array of height 2h, see Figs. |ljb) 
and (d). For large L this correspondence does not hold 
anymore, as demonstrated by comparison of panels (c) 
and (e) which display completely different symmetries. 
For dominoes, a mode with the same symmetry as the 
Yagi-Uda mode shown in panel (e) is forbidden because 
the horizontal electric field has to be zero at the metal- 
lic ground. This mode suppression is critical because 
the frequency of the fundamental mode of a Yagi-Uda 
structure, which is essentially controlled by the length of 
the long axis of the rods, is very sensitive to the lateral 
width L, being very different for short vertical rods and 
long horizontal rods. It is similarly important that the 
air gaps between the dominoes arc laterally open. Let us 
imagine that, instead of being open, they were laterally 
closed by PEC walls. Since the electric field inside the 
gap points longitudinally from one domino to the next 
one, it should then vanish at both PEC lateral walls, giv- 
ing a cosine- like x dependence. It is then clear that the 
mode frequency should increase when the lateral size L 
is decreased, arriving to the conclusion that the modes of 
such an structure have L-sensitive bands. In contrast, for 
DPs the boundary conditions at the laterally open gaps 
are more akin to perfect magnetic walls, and the corre- 
sponding field is approximately x-independent inside the 
air gaps even when L — > 0. 

Let us scrutinize in more detail the role played by the 
lateral dimension L, considering now realistic metals and 
paying attention to the spectral regime. The periodic- 
ity d is chosen to set the operating wavelength in the 
desired region of the EM spectrum, and L is varied in 
the range L = 0.5<i, . . . , 24c?, while the remaining param- 
eters are kept constant (a = 0.5c?, h — l.5d). Aluminum 
is selected for low frequencies, where metals behave al- 
most like PECs. In order to work at A = 1.6 mm, pro- 
viding an operating angular frequency of the order of 1 
THz, we first consider d = 200 /im. The evolution of the 
modal effective index as a function of the lateral dimen- 



sion normalized to the wavelength is plotted in Fig.JTJf). 
The curves corresponding to a PEC (black line) and alu- 
minum at A = 1.6 mm (red line) are, as expected, almost 
identical. We can now quantify the sensitivity of the ef- 
fective index to L, its variation being only about 12% 
even when L goes from L = oo to L = 0.5c? = A/16, well 
inside the sub- A regime. The rapid variation of n c g when 
L « d can be understood with the help of Eq. [T] In 
this limit the EM fields start to be z-dependent inside the 
air gaps, as verified by inspection of the fields (not shown 
here), letting us to conclude that q z ^ and consequently 
q y < ko- Then, if we assume that Eq. [Tl is still valid in 
this limit, it predicts that the modal effective index is 
reduced when L — > 0. To investigate the performance 
of DPs at higher frequencies, the structures have been 
scaled down by factors 1/10 and 1/100. The fact that 
the curves corresponding to A = 0.16 mm (green line) and 
A = 0.016 mm (blue line) do not lie on top of the previous 
ones is a signature of the departure of aluminum from the 
PEC behaviour. Nevertheless, even at A = 0.016 mm, the 
variation of the effective index is still smaller than 15%, 
a property that can be exploited for taper design in the 
far infrared, albeit with slightly reduced performance as 
compared with the THz case. When the operating fre- 
quency moves to the telecom regime (A = 1.5 /im, ma- 
genta line) the variation of the effective index is much 
larger (about 38%). Let us mention in passing that, for 
this wavelength, we have used gold because it has lower 
loss than aluminum, and that the geometric scaling had 
to be different (the periodicity is now d = 0.13 /im) be- 
cause with the same simple scaling as before no DP mode 
is supported due to the value of the permittivity of gold. 
The important message of panel (f) is that, although a 
variation of n e s begins to be noticeable when the lateral 
dimension L goes below A/2, DP bands arc fairly insen- 
sitive to L in the range L — 0.5c?, . . . , 24c? when operating 
at low frequencies. Such property has, to our knowledge, 
never been reported for geometrically-induced plasmon 
modes in corrugated wedges, channels, or wires. 

Now that we have identified the THz range as the most 
promising operation regime, we will characterize the cor- 
responding absorption (ohmic) and bending losses. Fig- 
ure|2ja) displays the propagation length [I = [2Im(fc)] _1 ) 
of DPs as a function of wavelength for various widths L 
(d = 200 /im for all curves). The only source of damp- 
ing in the domino structure is absorption by the metal 
as no radiation losses occur in such straight waveguide. 
Propagation length is very short for small A, when the 
Brillouin zone edge is reached, but it rises for increasing 
wavelength as the dispersion curve approaches the light 
line. Several lateral widths L are considered and it is ob- 
served that DPs in structures of smaller L travel longer 
distances before being damped. The radiation losses oc- 
curring in a waveguide of L = 0.5c? with a 90 degree 
bend are analyzed in panel (b). In this case the metal 
is modelled as a PEC for easier evaluation of radiation 
losses. Data for three radii of curvature of the bend, 
R = A/2,A,2A, are presented. In all three cases, losses 
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FIG. 2: Absorption (ohmic) and bend (radiation) 
losses, a, Normalized propagation length of DPs in rectilin- 
ear guides of various L as a function of A (ft = 1.5d, a = 0.5d, 
d = 200 /im). b, Bend loss of DPs for three radii of curvature 
as a function of A (h = 1.5d, a = 0.5d, L = 0.5d, d — 200 (im). 
Inset: Poynting vector field (modulus) distribution in a hor- 
izontal plane slightly above (30 /im) the height of the bend 
(top view). The chosen wavelength and radius of curvature 
are marked with an open black dot in panel (b). The red solid 
vertical lines in both panels indicate the operating wavelength 
used later, A = 1.6 mm. 



are low for small A but they grow as the wavelength is 
increased and the modal wave vector approaches fen. The 
dependence with the radius of curvature of the bend is 
not surprising either, radiation losses being higher for 
smaller radius of curvature. Consideration of absorption 
and bending losses suggests that A = 1.6 mm can be an 
appropriate operating wavelength. For A — 1.6 mm and 
L = 0.5d = 100 /im, the modal size is 0.34mm, i.e., about 
A/5. The propagation length is about 20A, much larger 
than the size of the devices presented later, and the bend 
loss for R — A is reasonably low, about 10%. Notice that 
for corrugated V-grooves, the bend loss for R = 2 A is 
around 50%G3. 

In the rest of the paper we present a number of de- 
vices enabling the spatial and temporal modulation of 
EM fields. The operating wavelength is A = 1.6 mm and 
d = 200 /im. In all ensuing computations the metal is 
considered as a PEC, to easily quantify t he device ra- 
diation loss. Waveguide taperingf 16 ! 18 ! 25 ! 2 ^ is interesting 
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FIG. 3: Subwavelength concentration of a domino 
plasmon. a, Poynting vector field (modulus) distribution 
in a horizontal plane slightly above (30 /im) the height of the 
tapered domino structure (top view). The lateral width is 
tapered from Li n = 16d to L out = 0.5d (h = 1.5d, a = 0.5d, 
d = 200 /im, A = 1.6 mm). b-e, Amplitude of electric field in 
transverse vertical planes (longitudinal views) at the locations 
shown by white dashed lines in a. The white bar in b showing 
the operating wavelength is valid for the last four panels. 



for field concentration and amplificat ion, b oth important 
requirements for imaging application d 27 * 28 * and circuit in- 
tegration^]. Tapering is expected to be easy in domino 
structures due to the insensitivity to lateral width L of 
the dispersion relation of DPs. This hypothesis is con- 
firmed in Fig. [3^a), which shows a top view of the prop- 
agation of power along a tapered domino structure. The 
lateral widths of the waveguide are L m = 16d at the input 
port (bottom) and L out = 0.5d at the output port (top). 
The length of the taper is 16d, i.e., two wavelengths at 
the operating frequency, corresponding to a taper semi- 
angle of 9 = 31°. Remarkably, reflection is smaller than 
2% and only 5% of the incoming power is lost as radia- 
tion due to the shrinkage of the lateral dimension along 
the taper. Panels (b)-(e) are cross sections at different 
positions along the taper, vividly showing the process of 
field subwavelength concentration and enhancement. As 
mentioned above, the same design does not work in the 
telecom regime due to the sensitivity of DP bands to the 
lateral dimension L at those frequencies. A more adi- 
abatic taper, i.e., with smaller change of n e ff per unit 
length, would not circumvent the problem because ab- 
sorption loss is high in the telecom regime (for instance, 
for A = 1.5 /im, d = 0.13 /im, and L = 16d, the propaga- 
tion length is I = 1.2 /im). 

We now discuss how t o construct power dividers and 
directional coupler^ 30 * 31 !. Figure [4] plots the dispersion 
relation of the modes supported by one and two paral- 
lel domino chains. Black solid line corresponds to the 
DP of L — 0.5d. For two parallel waveguides (both with 
L = 0.5d), interaction gives rise to symmetric and anti- 
symmetric supermodes which, for a separation s = 0.5d, 
are shown by magenta curves with the solid line depicting 
the even supermode whereas the dashed one represents 
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FIG. 4: Dispersion relations corresponding to one and 
two parallel domino structures. Black solid (blue dashed- 
dotted) line is for DP of L = 0.5d (L = 1.5d). Magenta 
solid (dashed) line is for the even (odd) supermode of paral- 
lel domino structures separated a distance s = 0.5d, whereas 
green solid (dashed) line is the corresponding supermode for 
s — 1.25d. The individual dominoes in double structures have 
L = 0.5d. The red solid horizontal line indicates the operat- 
ing wavelength, A = 1.6 mm, used later. Insets: longitudinal 
electric field for the odd and even supermodes of two parallel 
domino structures separated a distance s = 1.25d at the oper- 
ating wavelength, and displayed in a transverse cross section 
lying at the center of the inter-domino gaps. For all structures 
h = 1.5d, a = 0.5d, d = 200 ^m. 



the odd supermode. Blue dashed-dotted line corresponds 
to a wider structure of L = 1.5c? and, importantly, this 
band stays very close to that of the s = 0.5c? even super- 
mode. Due to this good overlapping, a simple linking of 
a wide domino structure with L = 1.5c? to two parallel 
domino chains (widths L = 0.5c? and separation s — 0.5c?) 
should not result in an important impedance mismatch. 
This is verified in Fig. [5^a), displaying a top view of a 
power divider based on the previous idea. The input port 
(L = 1.5c?) at the bottom feeds the wide domino which, 
after 5 periods, is linked to two parallel domino structures 
(L = 0.5c?, s = 0.5c?). Due to symmetry only the even 
supermode is excited. The lateral separation s between 
the individual waveguides is then increased, followed by 
90° waveguide bends. Simulations show that reflection 
is also negligible in this power splitter. Total radiation 
loss in the device is equal to that happening in the bends 
(previously discussed in Fig. |2jb)), each being 5% of the 
incoming power, which is small taking into account that 
the bend radius is equal to the operating wavelength. 

Coupling between parallel domino structures can also 
be easily understood with the help of Fig. |4j display- 
ing with green lines the even (solid) and odd (dashed) 
supermodes for parallel domino chains separated a dis- 
tance s = 1.25c?. The insets show the modal shape of 
both modes. A larger s is chosen in this simulation to 
prevent the cut-off of the odd supermode. The coupling 
length needed to exchange the carried power between the 
waveguides is estimated as 7r/|&cven — ^oddl — 2.05 mm. 
This value agrees with the result found in the simulation 



of a directional coupler. Figure [5jb) displays a top view 
of two parallel domino chains (L = 0.5c?) separated a dis- 
tance s — 1.25c?. A mode is launched in the left bottom 
waveguide and characteristic beatings with the predicted 
coupling length are observed when a second waveguide is 
located close to it. 

As a final example a waveguide ring resonator is 
demonstrated, able to route various frequencies to dif- 
ferent output portiP^. To show once more the robust- 
ness of our proposal, this device has been designed with 
posts of circular cross section as seen from above (ra- 
dius r — 0.28c?), instead of square cross section paral- 
lelepipeds. We have checked that similar results are ob- 
tained for the square geometry. The top view of this 
device is displayed in Fig.[5jc). For the chosen geometric 
parameters (ring radius R = 2 mm, minimum separation 
in the coupling sections s = 1.25c?), it is possible to drop 
a resonant wavelength A = 1.54 mm at the output port. 
Non-negligible radiation losses only occur at the bends in 
the waveguide ring and they are about 25% of the input 
power, reasonably small considering that no optimization 
was attempted. 

In summary, we have presented a new class of surface 
electromagnetic modes, termed domino plasmons, which 
feature an extraordinary characteristic: their insensibil- 
ity to the lateral dimension of the structure supporting 
them. We have shown how the guiding properties of these 
modes enable the planar routing of THz radiation at a 
deep subwavelength level. Finally, the flexibility and ver- 
satility of waveguides based on domino plasmons have 
been demonstrated through the implementation of a va- 
riety of functional devices, which may thus prove very 
useful as a new concept for subwavelength THz circuitry. 



III. METHODS 

All results in this paper have been obtained by means 
of numerical simulations performed with the Finite Ele- 
ment Method (FEM) using commercial (COMSOL Mul- 
tiphysics) software. A summary of the main modelling 
details follows. The technique is a volume FEM oper- 
ating in frequency domain. Due to the various length 
scales involved, highly non-uniform meshes are used. In 
order to ensure an adequate representation of the elec- 
tromagnetic fields, the size of the elements is a fraction 
of the skin depth inside the metallic parallelepipeds, a 
fraction of the characteristic geometric dimension (e.g., 
a) in the neighborhood of the structure, and a fraction of 
the operating wavelength at the boundaries of the sim- 
ulation domain. The mesh size is refined until the re- 
sults are stable with an accuracy of 1%. For instance, 
in the case of band structure computations, the conver- 
gence criterium is based on the value of both the real 
and imaginary parts of the wave vector. The final mesh 
is different for each simulation, but typical values of the 
tetrahedra sizes in different regions of the simulation do- 
main are of the order of 1/5 of the skin depth and char- 
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FIG. 5: Domino plasmon devices. Top view of: a, Power divider, b, Directional coupler, c, Waveguide ring resonator. The 
Poynting vector (modulus) is displayed in a horizontal plane slightly above (about 30 fim) the height of the domino structures. 
White bars show the operating wavelength in each panel. The various geometric parameters are described in the main text, 
the periodicity being d = 200 /im in all cases. 



acteristic geometric dimension, and 1/10 of the wave- 
length. The typical number of degrees of freedom lies 
between 3 x 10 5 and 3 x 10 6 , a range in which the ma- 
trix equations are inverted with direct solvers. For the 
computation of band structures, the corresponding eigen- 
value problem is posed in a single unit cell where Bloch 
boundary conditions together with Scattering boundary 
conditions are used. For device modelling, the excitation 
of the corresponding scattering problem is provided by 
Port boundary conditions, the ports being fed by the so- 
lutions of the corresponding eigenvalue problem. Open 
space is mimicked with Perfectly Matched Layers or/and 



Scattering boundary conditions. For the simulation of 
ideal metals, Perfect Electric Conductor (PEC) bound- 
ary conditions are employed. For realistic metals the 
dielectric permittivities are taken from Ref. [33J for Al 
(e A1 = -3.39 x 10 4 + «3.5 x 10 6 at A = 1.6mm), and 
Ref. [34] for Au (e Au = -103 + i8.7 at A = 1.5 //m). 
At low frequencies it was often sufficient to simulate the 
domains outside the conductors together with Surface 
Impedance Boundary Conditions, as verified by compar- 
ison with simulations fully accounting for the EM fields 
inside the metals. 
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